% This function runs the second step of the SR approach where macros appears
% linearly in the law of motion for the factors
function [theta2,AVarTheta2Step2] = step2SR_threshold_macro(AVarTheta1Step1,outputStep1,...
    bandwidthT,econAct,h0RegimeOn,hxRegimeOn,numBootStep2)

% option parameters
dispOn         = 1;
allowDeSaling  = 1;           % allow for descaling of the var_u to ensure omegeHat is positive definite
thresholdLevel = 0;
[theta2,deScalingX,deScalingZ] = Var1Threshold_GMM_closedForm_macro(...
    AVarTheta1Step1.VarFactorsRobust,...
    AVarTheta1Step1.autoCov10Factors,...
    AVarTheta1Step1.autoCov01Factors,...
    outputStep1.factors,outputStep1.macros',...
    econAct,thresholdLevel,allowDeSaling,h0RegimeOn,hxRegimeOn);
if dispOn == 1
    if deScalingX ~= 1
        disp(['The measurement error was descalted by ', num2str(deScalingX), ' in the estimated P dynamics for x'])
    end
    if deScalingZ ~= 1
        disp(['The measurement error was descalted by ', num2str(deScalingZ), ' in the estimated P dynamics for z'])
    end
end

% Ensure that the threshold model is stationary
[delta,theta2] = runInduceStatVAR_threshold_macro(theta2,...
    AVarTheta1Step1.VarFactorsRobust,...
    AVarTheta1Step1.autoCov10Factors,...
    AVarTheta1Step1.autoCov01Factors,...
    outputStep1.factors,outputStep1.macros',econAct,thresholdLevel,h0RegimeOn,hxRegimeOn);


if numBootStep2 > 0
    seedNum = 1;
    numBootStep1And3 = outputStep1.numBootStep1And3;
    disp('RUNNING BIAS-ADJUSTMENT OF P-DYNAMICS')
    outBiasCorrect  = biasCorrectTheta2_threshold_macro(theta2,AVarTheta1Step1.VarFactorsRobust,...
         AVarTheta1Step1.autoCov10Factors,...
        AVarTheta1Step1.autoCov01Factors,...
        outputStep1.factors,outputStep1.macros',econAct,thresholdLevel,...
        bandwidthT,h0RegimeOn,hxRegimeOn,allowDeSaling,numBootStep2,seedNum,numBootStep1And3);
    theta2 = outBiasCorrect.theta2Adj;
    delta  = outBiasCorrect.delta;
end

% The standard errors for theta2 using asymptotical inference
AVarTheta2Step2 = getAVarTheta2_threshold_macro(theta2,AVarTheta1Step1.VarFactorsRobust,...
    AVarTheta1Step1.autoCov10Factors,AVarTheta1Step1.autoCov01Factors,...
    outputStep1.factors,outputStep1.macros',bandwidthT,econAct,thresholdLevel,h0RegimeOn,hxRegimeOn);

% Additional output
AVarTheta2Step2.thresholdLevel = thresholdLevel;
AVarTheta2Step2.delta2         = delta;
end

